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Abstract 

The performance of collaborative beamforming is analyzed using the theory of random arrays. The 
statistical average and distribution of the beampattern of randomly generated phased arrays is derived 
in the framework of wireless ad hoc sensor networks. Each sensor node is assumed to have a single 
isotropic antenna and nodes in the cluster collaboratively transmit the signal such that the signal in 
the target direction is coherently added in the far-field region. It is shown that with N sensor nodes 
uniformly distributed over a disk, the directivity can approach N, provided that the nodes are located 
sparsely enough. The distribution of the maximum sidelobe peak is also studied. With the application 
to ad hoc networks in mind, two scenarios, closed-loop and open-loop, are considered. Associated with 
these scenarios, the effects of phase jitter and location estimation errors on the average beampattern are 
also analyzed. 

To Appear in IEEE Transactions on Signal Processing, 2005. 
I. Introduction 

Recent advances in the construction of low cost, low power, and mass produced micro sensors and 
Micro-Electro-Mechanical (MEM) systems have ushered in a new era in system design using distributed 
sensor networks [1,2]. The advent of sensor network technology provides a variety of applications that 
have been considered unrealistic in the past. One such application is in the area of space communications: 
with ad hoc sensor networks, a number of sensor nodes randomly placed on a planet can collaboratively 
collect information and then, also collaboratively, send the information back to Earth. In this scenario, 
the sensors must have an ability to transmit information over very long distances with high energy 
efficiency. In this kind of point-to-point communication, directional antennas are a preferred means to 
avoid interference. 

In general, this can be achieved by adaptive beamforming. Given a number of well-designed an- 
tenna elements at the transmitting/receiving sensor nodes, each node could in principle autonomously 
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in this publication are those of the authors and do not necessarily reflect the views of the National Science Foundation. 
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part at the IEEE International Conference on Acoustics, Speech, and Signal Processing (ICASSP'05), Philadelphia, PA, March 
2005. 
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Fig. 1. Collaborative beamforming concept in ad hoc sensor networks. 

transmit/receive the information to/from any desired direction. The advantages and applications of beam- 
forming with antenna arrays are well known; in wireless communications, this enables Space-Division 
Multiplex Access (SDMA), a technology which has the potential to significantly increase the capacity of 
multiple access channels. 

One of the most important constraints on wireless sensors is energy efficiency. Since the sensor 
nodes are often distributed in places where manual maintenance is costly, such as remote locations, 
on top of buildings and so on, it should be possible to operate these for several months without 
battery replacement. Considering the fact that each antenna element requires analog circuitry (and thus 
leads to costly hardware), in practice each distributed sensor is likely to be equipped with only a 
single antenna and this precludes the use of autonomous beamforming in scenarios of very energy 
efficient communication. Nevertheless, if the sensors in the cluster share the information a priori and 
synchronously transmit the data collaboratively as sketched in Fig.[2 it may be possible to beamform 
when transmitting (and also receiving) the data in a distributed manner. The resultant overhead due to 
intra-cluster information sharing may be relatively small as this can be done by low-cost short distance 
broadcasting-type communication among nodes. Thus, with distributed collaborative beamforming, the 
nodes can send the collected information to the far-end receiver over long distances with high efficiency. 
Also, only the sensor cluster in the specified target direction receives the data with high signal power 
and no significant interference occurs for clusters in other directions. Overall there is thus a potential 
to increase the capacity of the multiple access channel significantly despite the additional overhead for 
information sharing. 

The obvious question is whether one can form a nice beampattern with a narrow mainbeam, or 
achieve a reasonable directional gain (directivity). The sensor nodes in ad hoc networks are located 
randomly by nature, and the resultant beampattern depends on the particular realization of the sensor node 
locations. Therefore, it may be quite natural to treat the beampattern with probabilistic arguments. In this 
paper, assuming idealized channel model conditions and antenna properties, we analyze the achievable 



DRAFT 



3 



performance of collaborative beamforming based on distributed sensor nodes in a probabilistic sense. 
Specifically, the statistical properties of the achievable beampattern of the random sensor arrays are 
studied based on the following assumptions. The sensors form an ad hoc network and the geometry of 
the cluster is given by a two-dimensional disk of a given radius over which all sensor nodes are distributed 
uniformly as illustrated in Fig.^ Since the corresponding far-field beampattern depends on the particular 
realization of the random array of nodes, the probability distribution of the far-field beampattern is of 
particular interest. 

To the best of the authors' knowledge, the beampattern aspects of collaborative beamforming using 
random arrays have not been analyzed before in the framework of wireless ad hoc sensor networks. 
Nevertheless, in the antenna design literature, probabilistic analysis of random arrays is not new. In the 
framework of linear array design with a large number of sensors, Lo [3] has developed a comprehensive 
theory of random arrays in the late 1960s. It has been shown that randomly generated linear arrays with 
a large number of nodes can in fact form a good beampattern with high probability 1 and that with linear 
random arrays of N sensors, the directivity approaches N asymptotically. Although our scenario is quite 
different in that our main goal is not to design array geometry but to exploit the randomness of the 
distributed sensor network, it turns out that the results we shall develop in this paper can be seen as 
an extension of the theory of linear random arrays, [3], to random arrays on a disk. Thus, the same 
conclusion will be reached: with N collaborative sensor nodes, one can achieve a directivity of order N 
asymptotically. 

The major difference between classical beamforming by antenna arrays and distributed beamforming is 
that whereas the geometry of the former is usually known a priori, the exact location of the sensor nodes 
in ad hoc network is not, and it should be acquired dynamically. Even if their relative location is estimated 
by some adaptive algorithm (e.g., [1] for receiver beamforming), considering the low SNR operation of the 
sensor nodes, it is almost certain that the acquired geometric information has some inaccuracy. Also, since 
all nodes are operated with physically different local oscillators, each node may suffer from statistically 
independent phase offsets. In order to model and elucidate the effect of these impairments, we consider 
the following two scenarios: closed-loop and open-loop. The closed-loop scenario may be described as 
follows. Each node independently synchronizes itself to the beacon sent from the destination node (such 
as a base station) and adjusts its initial phase to it. Thus, the beam will be formed in the direction of 
arrival of the beacon. This kind of system is often referred to as a self-phasing array in the literature, 
and may be effective for systems operating in Time-Division Duplex (TDD) mode. The residual phase 
jitter due to synchronization and phase offset estimation among sensor nodes is then often the dominant 
impairment. On the other hand, in the open-loop scenario we assume that all nodes within the cluster 
acquire their relative locations from the beacon of a nearby reference point or cluster head. The beam will 
then be steered toward an arbitrary direction. Thus, the destination need not transmit a beacon, but each 
node requires knowledge of its relative position from a predetermined reference point within the cluster. 
This case may occur in ad hoc sensor networks where sensor nodes do not have sufficient knowledge 
of the destination direction a priori. In this scenario, since the acquisition of precise knowledge is not 
realistic, the effects of location estimation ambiguity among sensors upon the beampattern may be of 
particular interest. 

Throughout the paper, the nodes and channel are assumed to be static over the communication period, 
and for simplicity the information rate is assumed to be sufficiently low that Inter-Symbol Interference 
(ISI), due to residual timing offset, is negligible. It will also be assumed that all nodes share the same 
transmitting information a priori, as the main focus of the paper is on the beampattern, rather than the 

! It is interesting to note that the theory of random arrays has been discussed and developed almost exclusively in the antenna 
design community, e.g., in [3-6]. 
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Fig. 2. Definitions of notation. 

front-end communication performance. 

The paper is organized as follows. Section |H] describes the assumptions, model, and main parameters 
that describe beam characteristics associated with the framework of wireless ad hoc sensor networks. 
In Section the average properties of the beampattern are derived. The average beampattern of linear 
random arrays has been derived in [3], and our results can be seen as its extension to our sensor network 
model. For analytical purpose, we also introduce the concept of 3 dB sidelobe region. In Section IIVI 
the statistical distribution of the beampattern in a specific direction is derived. Lo [3] has derived the 
distribution of the beampattern in linear arrays based on a Gaussian approximation of the array factor, 
which is a common assumption in the random array literature. In contrast, we shall develop a numerical 
method to calculate the exact distribution of the beampattern and also examine the accuracy of the 
Gaussian approximation in detail. 

The distribution of the maximum of the sidelobe region is discussed in Section |V] This aspect of 
beampattern was analyzed by Steinberg [4], Agrawal and Lo [5], and Donvito and Kassam [6] in the 
framework of linear random arrays. In this paper, we derive an upper bound on the distribution of the 
maximum sidelobe in our framework of collaborative beamforming based on the approach of [6]. The 
effect of phase jitter or location estimation errors on the resultant beampattern, associated with the closed- 
loop and open-loop scenarios, is analyzed in Section Steinberg [7] has analyzed the effects of phase 
estimation errors in linear arrays, and based on a similar approach we analyze the effects of the average 
beampattern with phase estimation errors. Finally, Section IVlTl concludes the paper. 

II. System Model and Beampattern 

The geometrical configuration of the distributed nodes and destination (or target) is illustrated in Fig.|2] 
where, without loss of generality, all the collaborative sensor nodes are assumed to be located on the 
x-y plane. The kth node location is thus denoted in polar coordinates by (r^,^). The location of the 
destination is given in spherical coordinates by (A, <fto,9o). Following the standard notation in antenna 
theory [8], the angle 9 € [0, n] denotes the elevation direction, whereas the angle <fi G [— it, it] represents 
the azimuth direction. In order to simplify the analysis, the following assumptions are made: 

1) The location of each node is chosen randomly, following a uniform distribution within a disk of 
radius R. 

2) Each node is equipped with a single ideal isotropic antenna. 

3) All sensor nodes transmit identical energies, and the path losses of all nodes are also identical. 
Thus the underlying model falls within the framework of phased arrays. 
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4) There is no reflection or scattering of the signal. Thus, there is no multipath fading or shadowing. 

5) The nodes are sufficiently separated that any mutual coupling effects [8] among the antennas of 
different sensor nodes are negligible. 

Furthermore, we also assume that all the nodes are perfectly synchronized so that no frequency offset 
or phase jitter occurs. The effects of phase ambiguities on the beampattern will be discussed in Section lVll 

Let dk(4>,0) denote the Euclidean distance between the kth node and the reference location (A, (/),$), 
which is written as 



d k (0,6) = ^A 2 + rl-2r k Asm6cos((i>-?p k ) . (1) 
If the initial phase of node k € {1, 2, . . . , N} is set to 

^ = -y4(^o,fio), (2) 

the corresponding array factor, given the realization of node locations r = [n, r<i, ■ ■ ■ , r^] G [0, R] N and 

ip = [ipi,i>2, ■ ■ ■ , V'iV"] ^ [ -7r > 7r ] Af > is written as 

N N 

F{M\t,iI>) = i^ e ^f<W) = I£ e *?[*(W-i(*..ft»)] j ( 3 ) 

k=l k=l 

where N is the number of sensor nodes and A is the wavelength of the radio frequency (RF) carrier. 

In this paper, we are interested in the radiation pattern in the far-field region, and we assume that the 
far-field condition A^> r k holds. The far-field distance d k {4>, 0) in Q can then be approximated as 

d k (4>,0) « A - r k s'm6cos((j) - ifj k ). (4) 

The far-field beampattern is thus approximated by 



F(<j>,e\r,iJ)) « l^ e i¥^[ sin9 ocos(0 o -^)-sinecos(^ fc )] ± p ((j> , 9\r , . (5) 
fe=i 

Alternatively, instead of applying ^ k as in if we choose 



+ 2ir 

® k = ~^r k sin #0 cos(0 o - W, (6) 



then we obtain the array factor as 

N 



F\<j),e\r,iP) = 1 JV*ie^(^) 



k=l 
N 



_ ST^ f>3 2 ri A ~ rk sin6»cos(0-t/) t! )+r fc sin0 o cos(^ — 4>k)} 



N 
fc=i 



j2z-a}_ e j^rt.[sinfl eos(<ft — 0»)-sinflcoB(4>— 4> k )] A (<fr, 6\r,tf>). 



(7) 



k=l 



The only difference between P((j),6\r,rp) in (I5} and F*(<f>,9\r,iJ)) in Q is the existence of the initial 
phase offset of ^r-A The far-field beampattern is thus identical for both systems, and the received signal 
exhibits no difference as long as the base station compensates for this phase rotation. 

Therefore, there are two ways of forming a beam. One way is to use (0), but this approach requires 
accurate knowledge of the distance, relative to the wavelength A, between each node and the destination. 
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Thus, this applies to the closed-loop case such as self-phasing arrays. Alternatively, the use of © requires 
knowledge of the node positions relative to some common reference (such as the origin in this example), 
and thus corresponds to the open-loop case. Knowledge of the elevation direction 9q is also required, 
but this may be assumed to be known a priori in many applications. In both cases, the synchronization 
among sensors is critical, which may be achieved by the use of reference signals such as those of the 
Global Positioning System (GPS). 

Of particular interest in practice is the case where 9q = |, i.e., the destination node is in the same 
plane as the collaborative sensor nodes. Therefore, we will consider the beampattern in this plane and 
thus assume that 9 = 9q = | for the rest of the paper. We then denote F((f>,9 = ir/2\r,ifr) in © by 
F(4>\r,tf>) and F^(4>,9 = 7r/2|r, if>) in Q by (4>\r , ?/>) for simplicity. 

By assumption, the node locations (r k , ipf.) follow the uniform distribution over the disk of radius R. 
Thus, the probability density functions (pdfs) of r k and ip k are given by 

2r 1 
fr k ( r ) = ft2> 0<r<i?, and fy^) = — -vr < ip < vr. 

From ©, we have (with 9 = 9o = 7r/2) 

N N 
k=l k=l 

where f k = r k /R and tjjk — ipk — The compound random variable 

z k = f k s'mip k , (9) 

has the following pdf: 

f Zk (z) = -Vl-z 2 , -1<z<1. (10) 

7T 

Note that since the above model is symmetric with respect to the azimuth direction cj), any particular 
choice of 4>q does not change the results in the following. Therefore, without loss of generality, we 
assume that <po = 0, and the parameter <fi simply corresponds to the difference angle between the target 
direction and the reference. Also, note that \(p\ <ir. 
The array factor of (JSJl can then be rewritten as 

1 N 

F(0\ z ) = jjY, e ~ jAwRsin{i)zk > < n ) 

k=l 

where R = y is the radius of the disk normalized by the wavelength. 
Finally, the far-field beampattern can be defined as 

P(<t>\z)±\F(<t>\z)\ 2 = F{<p\z)F*{4>\z) 

N N 

= _l ^ e~ j47rRs ' m ^y zk ~ z,) 
^ k=i i=i 

N N 

k=l 1=1 
l^k 



where 



a(<t>) =47ri?sin|. (13) 
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Fig. 3. Average beampattem with different R and N — 16, 256. 
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III. Average Properties of Beampattern of Uniformly Distributed Sensor Array with 

Perfect Phase Information 

A. Average Far-Field Beampattern 

We start by investigating the average beampattern of the random array resulting from the distributed 
sensor network model in the previous section. Here, the average is taken over all realizations of z, and 
from (112ft the average beampattern is expressed as 

P m (<l>)±E g {P(<f>\z)}, (14) 

where E x {-} denotes expectation with respect to the random variables x. From ill 2b and (flOl l. it can be 
readily shown that 

Ji W)) 

a(4>) 

where J n (x) is the nth order Bessel function of the first kind. Although the function Ji {x)/x is oscillatory, 
the local maxima of oscillation tend to decrease with increasing x. In (I15t . the first term represents the 
average power level of the sidelobe, which does not depend on the node location, whereas the second 
term is the contribution of the mainlobe factor. Since, conditioned on tp, the array factor of the form (fTTT i 
is an average of bounded independent and identically distributed (i.i.d.) complex random variables, by the 
weak law of large numbers the beampattern (fT2l converges to the ensemble average (fT3l in probability 
as N — > do. 

The average beampattern of dl 5b is plotted in Fig.|3]for several values of R with N = 16 and 256. As 
can be observed, the sidelobe approaches 1/N as the beam angle moves away from the target direction. 
To gain further insight, consider the asymptotic expansion of the Bessel function J\ (x) for x S> 1 as 

JAx) ~ \ COS { X -7T ) . (16) 

V 7TX \ 4 / 



DRAFT 



s 



We then have 



Ji(x) 



■ cos 



TTX- 



-IT 



and (fl5b becomes, for a{4>) = AttR sin ( | ) » 1, 



cos 



«(0) - T 71 " • 



The nth peak of the average sidelobe will appear around a(4> n ) as (n — ^) 7r, 
corresponding value becomes 



n 



1,2, 



7r n 



i) 



(17) 

(18) 
and its 

(19) 



which does not depend on R. The nth peak and nth zero positions (in the sense of the second term in 
(fT5l ) can then be expressed asymptotically as 



</>P eak ~2arcsin 



n 



2arcsm ^ 

\ AR 



(20) 



(21) 



Since the peak sidelobe value does not depend on R and is less sensitive to the value of N, it is apparent 
that the only way one can avoid high peaks in the sidelobe region is to increase R such that most of the 
major peaks are relatively concentrated around the mainlobe. This phenomenon will be further examined 
in the following subsections. 

B. 3 dB Beamwidth of the Average Beampattern 

One of the important figures of merit in directional antenna design is the 3 dB beamwidth. In the 
deterministic antenna, the 3 dB beamwidth is the threshold angle at which the power of the beampattern 
drops 3dB below that in the target direction <pQ. In our scenario, the 3dB beamwidth itself is a random 
variable and it is not easy to characterize analytically. Thus, as an alternative measure, we may define 
the 3dB beamwidth of the average beampattern denoted by as the angle <j> that satisfies 

I- (22) 



t3dB^ 



In the limit as N — > oo, one may obtain 



i 3dB 



2 arcsin 



0.1286\ 
R J 



(23) 



by numerically solving (1151 . For R 3> 1, d23l > can be approximated as (j)f^ B pa 0.26/-R. Therefore, the 
beamwidth is asymptotically independent of N and is mainly determined by the inverse of the disk radius 
of the cluster. Consequently, sparsely distributed sensors form a narrow beam on average provided that 
the cluster radius is sufficiently large. 

This sharp mainbeam property may be desirable, but if the far-field destination node has mobility, it 
should be designed carefully; the calibration should take place before the mobile node moves out of the 
mainbeam, but the mainbeam width is inversely proportional to the normalized radius R as observed in 
Fig-El Therefore, calibration should be performed more frequently if the destination node moves rapidly 
or when R is increased. 
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Fig. 4. Definitions of 3 dB beamwidth and 3 dB sidelobe region. 



C. 3 dB Sidelobe Region 

Similar to the 3 dB beamwidth concept, it may also be convenient for our subsequent analysis to define 
the region within which the average of the sidelobe beampattern falls below some threshold level. As 
we have seen, for large R, the sidelobe of the average beampattern is given by 1/N asymptotically. 
Therefore, we shall define the 3 dB sidelobe region as the region in which neither neighboring sidelobe 
peak in the average beampattern exceeds 3 dB above 1/N. Let no denote the minimum index of the peak 
position such that the corresponding peak value satisfies this 3 dB condition. Specifically, from H91 . no 
is the minimum integer n that satisfies 



,peakx 



1/N 



and it can be bounded by 



~1 + (JV-1) 



1 2 

n > - + - 

4 7T 



7T 



<2, 



N - 1 



7T 



(24) 



(25) 



Let 



'n 



> denote the angle corresponding to the zero point next to the noth peak sidelobe which can 



be obtained from d2lT) with 



n 



n . Consequently, in this paper, the 3 dB sidelobe region is defined as 

5 3 dB = {0|vr>|<A|>C'°}- (26) 

Fig.|H illustrates the definitions of the 3 dB sidelobe region together with that of the 3 dB beamwidth. 

As will be shown in Section HV-CI the idea behind the introduction of 3dB sidelobes is that in this 
region one may assume that the mean value of the random array factor of © sampled at cj> G 
becomes a random variable with approximately zero mean. Thus we may reasonably assume that the 
array factor has zero mean in this region, and this assumption significantly simplifies the analysis of the 
statistical distribution in the following sections. 
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Fig. 5. Threshold of 3dB beamwidth and 3dB sidelobe region with respect to R and number of nodes N. 



Fig-El shows the threshold angle above which the 3 dB sidelobe region begins. The asymptotic 3 dB 
beamwidth d23b is also shown for reference. As can be observed, whereas the 3dB beamwidth is less 
sensitive to the number of nodes N, the 3 dB sidelobe region will be considerably reduced as N increases. 
This means that as N increases the dominant non-negligible sidelobe peak may occur with high probability 
unless R is also increased. This trade-off will be clarified by the study of directivity in the following 
subsection. 



D. Average Directivity 

The directivity or directional gain is the parameter that characterizes how much radiated energy is 
concentrated in the desired direction, relative to a single isotropic antenna. Specifically, it may be defined 
as 



D 



r p 

J—-K 



r p 

J —TV 



(27) 



where P(<p) is the radiated power density in the direction of (f>. In the scenario of this paper, since P{4>) 
depends on the particular realization of z, the corresponding gain may be expressed, by substituting 
P((p\z) of (11 2b into the above, as 



D(z) 



1 1 N N 



fc=i i=i 



The mean value of (128b is given by 



L> av 4 E z {D{z)} 



(28) 



(29) 
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Fig. 6. The relationship between directivity D m /N and normalized radius _R. 
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Unfortunately, direct calculation of d29t does not result in a closed-form or insightful expression. Thus, 
we shall consider the following as an alternative measure: 

2vr 2vr 



Note that by Jensen's inequality, we have 



(30) 



(31) 



which means that D m in (l3*0l is a lower bound on D m . However, since by the law of large numbers the 
denominator of D(z) may approach its average value with high probability as N increases, the above 
bound is expected to become tight as N increases. This will be verified in the numerical results below. 
Substituting dT3l into d30l . we obtain 

(32) 



1 + (N - 1) 2^3 



1 3.i no. 

2 j 2 ' ' ' ' 



-(4vri?) : 



where 2-F3 (|, § ; 1, 2, 3 ; — x 2 ) is a generalized hypergeometric function which monotonically decreases 
with increasing x and converges to as x — ► 00. Therefore, unlike well-designed deterministic linear 
arrays, the gain of a given realization is very likely to be less than N, and the limit N can be approached 
only by increasing R. This agrees with the previous observation that the average mainbeam becomes 
narrow as R increases and thus improves the directivity. 

It should be noted that although (l3*2l has a simple form and offers some insight into the asymptotic be- 
havior of directivity, the calculation of the generalized hypergeometric function involved in d32t becomes 
numerically unstable as R increases, and it is much easier to numerically integrate the denominator of 
d30l directly. 
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Fig. 7. The relationship between directivity D w /N and node density N/R. 



Fig-El shows the relationship between the normalized directivity bound D w /N and R. Also shown in 
the figure as diamond-shaped points are the corresponding exact average directivities D av /N obtained by 
the simulation of 1000 realizations. As can be observed, the bound is very tight compared to the exact 
performance. Thus, it follows that in order to achieve high normalized directivity (i.e. directivity close 
to N) with N nodes, the distribution of the nodes should be as sparse as possible. In fact, we have the 
following theorem: 

Theorem 1 (Normalized Directivity Lower Bound): For large R and N, D av /N is lower bounded by 

£*L > (33) 

where \i is a positive constant independent of N and R (/x « 0.09332). 

Proof: See Appendix U ■ 

Note that the factor N/R which appears in d33l can be seen as a one -dimensional node density. To 
verify the above theorem, Fig.0 shows the relationship between D av /N and the node density N/R, as 
well as the lower bound in <33b . 

The above theorem indicates that there is a simple relationship between directivity and node density. 
It can be seen that the node density almost uniquely determines the normalized directivity D ay /N. It is 
important to note that in order to achieve a certain normalized directivity with large numbers of nodes 
N, the node density should be maintained to the desired value by spreading the nodes as sparsely as 
possible. Alternatively, if the normalized region R is fixed, it is not efficient in terms of normalized 
directivity to increase the number of sensor nodes. 

The above theorem also indicates that if the sensor nodes are uniformly distributed and if we are to 
choose N nodes out of them, in terms of normalized directivity it may be better to choose them as 
sparsely as possible. 
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Fig. 8. Average and realization of beampattem with R = 2 and N = 16. 

IV. Distribution of Far-Field Beampattern of Collaborative Beamforming with 



In the previous section, we have seen that random arrays have nice average beampatterns with low 
sidelobes. However, the average behavior does not necessarily approximate a beampattern of any given 
realization unless N — > oo. In fact, even though the average beampattern has a sharp mainbeam and 
sidelobes always close to l/N, there is a large dynamic range of sidelobes among randomly generated 
beampatterns. As an example, the average beampattern and one particular realization of randomly gen- 
erated sensor locations is shown in Fig.|8] The mainbeam of the realization closely matches the average, 
but sidelobes may fluctuate with a large dynamic range and easily exceed the average level. 

Therefore, in practice, the statistical distribution of beampatterns and sidelobes in particular, is of 
interest. By approximating the beampattern sidelobes as a complex Gaussian process, Lo [3] has derived 
the distribution of the beampattern in the case of linear random arrays. 

In the following, we first derive a numerical method that allows calculation of the exact distribution 
of the beampattern without applying Gaussian approximations. We then derive a convenient asymptotic 
form of the sidelobe distribution using a Gaussian approximation similar to [3], and evaluate its validity 
in our framework. 

A. Exact Evaluation of Distribution 

Since the array factor is a sum of i.i.d. random variables, its distribution can be computed numerically 
by the characteristic function method. To this end, from (11 Q let 



Perfect Phase Information 




(34) 



k=l 
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where 



x k = cos {z k a{4>)) , y k = sin {z k a(<j))) 
and a{4>) is defined in dT^b . The joint characteristic function of x k and y k is written as 

j[u cos(z k a(4)))+u sin(« fc a(^))] 1 



(35) 



(36) 



For a given pair of u and i/, the above expectation is a single integral of a well-behaved function and 
can be calculated numerically. 

Since F(cj)\z) is a sum of N i.i.d. complex random variables, the joint probability density of X and 
Y in d34b is given by 

(i \ 2 poo />oo 

The above integral can be computed efficiently using the two-dimensional Fast Fourier Transform (FFT). 
Finally, the complementary cumulative distribution function (CCDF) of the beampattern, i.e., the proba- 
bility that the instantaneous power of a given realization in the direction exceeds a threshold power, 
Pq, is given by 

^X 2 + Y 2 



Pr [P(<P) > P ] = Pr 



N 2 



>Pq 



fxv( x >y) dx d y- 



(38) 



x 2 +y 2 >N 2 P 



B. Gaussian Approximation of Distribution 

The exact evaluation of the CCDF outlined above is computationally demanding, especially if the 
desired numerical precision is high. Considering that the array factor consists of a sum of N statistically 
independent random variables, as N increases, by the central limit theorem we may expect that the array 
factor with any given direction, except at the deterministic angle = 0, approaches a complex Gaussian 
distribution. This approximation may typically result in a simpler distribution formula. To this end, we 
write d3*4l as 

F(<f>\z) = ^=(X - jY) (39) 



N 



where 



X 



1 - 

— E 



N 



cos \z k a 



M)' ^ = ^I>in(z fc a(0)). 
VN k=i 



(40) 



Since the z k 's are i.i.d. random variables, as iV increases the distribution of X and Y at the direction 
7r > |0| > may approach that of a Gaussian random variable with 



E{X} 



Var(X) 



t Ji W)) 

at(4>) 



N = TTLrr 



1 + 



E{Y} = 

1 

2 



Var(y) = - 1 



Ji (2«(0)) 
«(0) 

Ji(2a(0)) \ 
a{4>) ) 



W )) 

a(0) 



n 2 



A ^-2 



A ~ 2 
= V 



(41) 

(42) 
(43) 
(44) 
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Note that E{XY} = 0, i.e., X and Y are orthogonal and thus statistically uncorrelated. The joint pdf 
of X and Y is then given by 



fxy{x,y) 



1 



2na x a y 



exp 



\x — m a 
2al 



(45) 



The CCDF of Pq can be expressed as 

-X 2 + Y 2 



Pr [P(if>) > P ] = Pr 



N 



>Po 



Pr 



T 1 v 2 



exp 



Vx 2 + Y 2 > ^/WF 

2 



r coso; — m. 



2a 2 



9 • 2 

r sin a; 



da; dr 



where 



A-Ka x a y U 2 



I cos 2 a; sin 2 u> 



y/nV u erfc(W u -V u ) + e 



-(W^-K,) 2 



2a 2 



+ 



m x cos a; 

2a 2 L/ aJ 



dw, (46) 



(47) 



For a(0) S> 1, the terms Ji(2a((j>))/a(cj>) and | Ji(a(0))/a(</>)| 2 in the variance expressions d42l and 
d4*4*l) rapidly decrease and their contribution to the resulting variances becomes minor. Therefore, it is 
very likely that both variances are approximately equal in the sidelobe region. When this is the case, 
i.e., if cr 2 « cr 2 1/2, the distribution of the complex envelope becomes a Nakagami-Rice distribution. 
Consequently, the resulting integral can be expressed in terms of the first-order Marcum-Q function 



Pr [P(0 >P ] = Q = Q (V2m x , V2~NP 



(48) 



Furthermore, if the mean E {X} is zero, the envelope follows a Rayleigh distribution and we simply 
have 



Pr [P(cP) > P ] = e 



-NPn 



(49) 



C. Mean Value of Array Factor within 3 dB Sidelobe Region 

As we have seen, if the mean value of the array factor can be assumed to be zero, the distribution 
can be significantly simplified and thus analysis becomes readily tractable. From d4Tb it is apparent that 
under the constant variance constraint the mean value increases as N increases. Therefore, when N is 
large, the zero mean assumption may not be guaranteed in general. In Section ITlI-CI we have introduced 
the 3 dB sidelobe region, and in the following we derive properties of the mean value of the array factor 
in this region. 

From the definition of (1246 . the sidelobe in the 3 dB region satisfies 

NPM < 2. (50) 

It follows that 

Var(X) + Var(y) + \E {X}\ 2 < 2. (51) 
From (EJ, d42l . and d44t . we have 

1-1\E{X}\ 2 + \E{X}\ 2 <2, (52) 
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and thus we get 

2 



E{X}\ Z <- -. (53) 



1-F 



Therefore, in the 3 dB sidelobe region, the square of the mean is bounded by unity in the large-iV 
asymptote and thus the mean does not grow unbounded with the number of nodes N. 



D. Numerical Comparison 

In Fig.|9l the CCDF's computed with various formulae are shown with R = 2 and <j> = 7r/4, which 
corresponds to the sidelobe region. The exact formula of d3*Ht . the equal variance Gaussian approximation 
of d48l . and the zero-mean Gaussian approximation of d49l are shown in the figure. Note that the precise 
Gaussian case of d46l was also calculated by numerical integration, but it is almost identical to the 
Marcum-Q function approximation in d48t for this case and thus is not shown. As observed from Fig.|9l 
even the zero-mean Gaussian approximation may be valid for this sidelobe region, but for N = 1024 the 
Gaussian approximation will have some noticeable discrepancy with the exact value. This is due to the 
fact that the zero-mean approximation does not hold for this case. In fact, Fig.|5] indicates that for this 
value of N, the angle falls between the 3 dB sidelobe region and the mainlobe region and thus the zero 
mean assumption may not be accurate. 

Fig.EH snows the distribution at 3dB beamwidth of the average beampattern defined by d2"3l . In this 
case, the exact form d38l ». the precise Gaussian d46l . and the Marcum-Q approximation (I48t show different 
results even for relatively large N, since at this angle the variance of the array factor is small and a large 
number of random variables must be summed in d40l for a non-zero mean Gaussian random variable to 
adequately approximate the sum. As observed, as the number of nodes increases, the mainbeam variance 
becomes small and approaches the mean value of -3 dB, as expected. Therefore, for large N the mainbeam 
can be made stable. This observation agrees with a similar result for linear arrays in [3]. 
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Fig. 10. Distribution of beampattern at <j> = (j>^ with exact, precise Gaussian, and Marcum-Q formulas. 



V. Distribution of the Maximum of the Sidelobes 

It is well known that unlike periodic or equally-spaced antenna arrays, many arrays with unequal 
spacing will yield no grating lobes for a large number of elements. This property is also preserved for 
random arrays [3,7], but in order to verify this, one may need to find the distribution of the maximum 
power of the sidelobes. In this section, we develop an approximate upper bound on the distribution of 
the peak sidelobes for random sensor networks. 

In the previous section, we have seen that the distribution of the beampattern samples within the 
sidelobe region can be characterized by a zero mean Gaussian random variable if the zero-mean condition 
is satisfied. In the following, we further assume that the beampattern is a Gaussian random process. In this 
case, any two samples taken from the beampattern should be characterized by jointly Gaussian random 
variables. In the linear random array framework, the distribution of peak sidelobes has been studied in 
[4-6], assuming the array factor is a Gaussian process. For simplicity, only the 3dB sidelobe region is 
considered and it will be assumed that the process is stationary with zero mean. The extension to the 
non-stationary case is studied in [6]. 

In the following, the CCDF of the maximum peak sidelobe, which is the probability that the maximum 
peak sidelobe exceeds a given power level, will also be referred to as outage probability and denoted by 

-fout- 



A. Upper Bound on the Distribution of Peak Sidelobe 

Let v(a) denote the random variable representing the number of upward crossings at a given level a 
per interval in the 3dB sidelobe region 53dB- As shown in Appendix HT1 assuming that the array factor 
in this region can be approximated as a zero-mean Gaussian process, the mean of u(a) is given by 

E{y{a)} = 4 (l - sin J ^Rae~ a \ (54) 
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Fig. 11. Comparison of CCDF and upper bound of the sidelobe peaks with the node density N/R = 2. 



Note that the above function monotonically decreases with increasing a only for a > 1/ v2, and thus is 
meaningful only in this region. Finally, noticing that the outage probability is the probability that at least 
one peak exceeds level a and is given by [6] 

oo oo 

P out = Pr [v(a) > 1] = Pr[i/(a) = k] < ^ k Pr[v(a) = k) = E {v{a)} , (55) 

k=l k=l 

then d54t serves as an upper bound for the outage probability for the maximum sidelobe peak for 
a > 1/V2. Thus, we obtain the CCDF upper bound as 

1 



Pr 



maxl 2 + r 2 > P 



< 4 1 - sin 



izero 
J n 



NPn 



for NP > 



(56) 



B. Numerical Results 

Fig.^2 shows a comparison between simulation results and the upper bound (1561) . For the simulation, 
the outage probability is calculated based on 10000 randomly generated realizations with the node density 
N/R = 2 and only the peaks within the 3dB sidelobe region are examined. Also, in order to capture 
peak values accurately, the entire 3dB sidelobe region of (ft is sampled at a rate as large as WttR. As 
can be observed, the bound is in good agreement with simulation for large N. 

Let Po = NPq denote the threshold of the maximum peak value (Pn) normalized by the average 
sidelobe level (1/N). Since from (12 1 b <j?" o ° approaches zero as R increases, d56T l reduces to 



Pout < V Poe- Po , Po > 1/2. 



(57) 



The above inequality illuminates the relationship between the outage probability and R (assuming that 
(ftn™ is negligibly small). Fig.^] shows the maximum possible value of Po for a given outage probability 
and R. As can be observed, the maximum sidelobe may grow as R increases, but the amount is below 
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Fig. 12. Bound on sidelobe maximum with a given outage P out . 



12dB for many cases of interest. Consequently, the maximum sidelobe level (in the 3dB region) of 
randomly generated arrays may be written as Pq = Pq/N, where the required margin Pq depends on 
the parameters R and P out , but not on N. Thus, increasing N always results in a reduction of maximum 
sidelobe level. 



VI. Performance of Distributed Beamforming with Imperfect Phase 

So far, we have evaluated the beampattern assuming perfect knowledge of the initial phase for each 
node. In this section, we analyze the effect of the phase ambiguities in the closed-loop scenario as well 
as location estimation errors in the open-loop scenario. For each of the two scenarios, we derive the 
average beampattern and calculate the amount of mainbeam degradation. 



A. Closed-loop Case 

In the closed-loop case, the effects of imperfect phase may be easily derived, following the approach 
developed by Steinberg [7]. The initial phase of node k in (|2j will now be given by 

** = -ydfc(0o,0o) + ¥>* (58) 

where tp^ corresponds to the phase offset due to the phase ambiguity caused by carrier phase jitter or 
offset between the transmitter and receiver nodes. In the following, the phase offset tp^s are assumed to 
be i.i.d. random variables. Then, from ©, ©, and il 11 1, the far-field array factor (with 9 = 9q = tt/2) 
will be given by 

TV N 

F((f>\z,cp) = — ^ e^ - ** 4 **™ !+*".) = _J2e- jZk ^ Rain U j ' p \ (59) 

k=l k=l 
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The average beampattern of dT4b will be replaced by 

P a M = E z , lp {P(<l>\z,cp)}. 
Similar to (fT5l . direct calculation of d60t results in 



J! {a{4>)) 



a(4>) 



where 



(60) 



(61) 



(62) 



Thus, as N —> oo, the average beampattern will simply become a version of the original scaled by a 
factor of \A ip \ 2 . 

Let us now assume that the phase offset follows a Tikhonov distribution, a typical phase jitter model 
for phase-locked loop (PLL) circuits given in [9], 

1 



exp (cos(x)/a ip ) 



\X\ < IT, 



(63) 



2^o (1/4) 

where is the variance of the phase noise and I n is the nth order modified Bessel function of the first 
kind. The corresponding attenuation factor is given by 

v WHY 

The variance of the phase noise cr 2 is related to the loop SNR of the PLL by 



(64) 



(65) 



Fig. 1131 shows the degradation factor \A ip \ 2 with respect to the loop SNR. As observed from the figure, 
a loop SNR of at least 3 dB may be necessary for each node in order to reduce the overall beampattern 
degradation to less than 3dB. 



B. Open-loop Case 

In the open-loop case, our model of the initial phase is given in © with 9q = § > and if there are 
estimation errors in the location parameters r k and ipk, the initial phase will be replaced by 



$1 = y ( r fc + Sr k) cosOo - (tpk + Sipk)) 

—r k cos((/>o - (ipk + Sipk)) + —$ r k cos(0 o - {tpk + Sipk)), 



2ir 



(66) 



where 5r k and 5ipk are the corresponding error random variables, each set assumed to be i.i.d. and also 
independent of r k and ipk for simplicity. With the far-field approximation, we have 

y dk (<P, |) + * I ~ y {^4 - r k [cos( 



ipk) ~ cos(0 o -ipk~ Sipk)] + Sr k cos(</> - (ipk + Sipk))} 



— A + — r fc sin I Vfc 

+ y^fe cos(V>fe - Oo - 



8111 



Siph 



(67) 
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Fig. 13. Mainbeam degradation due to the phase noise in the closed-loop scenario. 



Let ip k ^ijj k - i+is^h , Then, the right-hand side (RHS) of © is given by 

2vr 4vr . ~ . / 4> - cf) - 5ip k \ , 2vr / ~ 6 - (j) + 5ip k 
— A - —r k sin tp k sin I - I + -y <^ fc cos lip k -\ - 

The resulting far-field array factor of Q will then be given by 

r N 

F(<j>\r, i/,, Srf,, Sr) = e^ A - £ sin ^ ™(^F**)+i^r» 

k=l 

and the beampattern is expressed as 



N N 



p(<p\z,v,s^) = - + i.rr e -^'K 1 ^)-*-»( tl ¥ iil )}e»-? 



(Vk-Vl) 



k=l 1=1 



where 



' k s'mtp k = f fc sin ( ip k + 



R 



v k = 5r k cos 4> k + 



2 2 



5r k cos (ip k + 5ip k - 4> ) . 



(68) 



(69) 



(70) 



(71) 
(72) 



Conditioned on 0, 0o and 5tp k , the angle ip k can be seen as a uniformly distributed random variable, and 
thus the pdf of z k is given by d!0l >. Considering the fact that r k and 5r k are assumed to be statistically 
independent, we further assume for analytical purposes that z k and v k are statistically independent. Then, 
again, the beampattern does not depend on the particular choice of <pQ. Furthermore, on modeling 5r k 
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as being uniformly distributed over [— r max , r max ] and assuming the phase term of Vk to be uniformly 
distributed over [0,2tt], the probability density function of Vk will be given by 



vrr n 



\v\ < r max . (73) 



Consequently, the average beampattern can be written as 



PM = ^ + (i - ^ ) \M<t>)\ 2 IAI 2 , (74) 



where 



A T = E n [e'f' ! - } = IJ cos (y w) In 



1 + 



- ,j, ,„ i . | Ji 1 4-irR sin *^~^ fc 

* {^-K^)} _ ^ j , <V6) 

and without loss of generality <po = was assumed. In d75t . \Fz (5 ; 1, | ; — x 2 ) denotes a generalized 
hypergeometric function which has an oscillatory tail but converges to zero as x increases. 

Also, assuming that the dipk are uniformly distributed over [— ^ max , VWx] and using the approximation 
sin (cj) + Sipk ) ~ <fi + dipk which is valid for the beampattern around the mainbeam, we obtain 

M<f>) ^U 1 ' ^) ^ Q; §,2; -(7T^(0 + Vmax)) 2 

1 A M „ /l 3 . , , , 2 



+ H 1+ ^J lFa U ; 2' 2; -^-^J- (77) 

Since the hypergeometric function 1F2 (| ; |, 2 ; — x 2 ) has a maximum peak value of 1 at x = 0, the 
above expression indicates that regardless of the value of R, there may be two symmetric peaks around 
the mainbeam at <f> = ±?/; max resulting in a pointing error. Therefore, the mainbeam may spread over by 
a factor of VWx- At the center of the mainbeam, we have 

^(0)=iF 2 ^;|2;-(vr^i)^. (78) 

Fig. [T?l shows the degradation factor \A r \ 2 and |^(0)| 2 for a given ^ and '^ ma * . As observed 
from the figure and discussion above, the angle estimation error has two effects, i.e., pointing error and 
mainbeam degradation. In particular, if we wish to suppress the mainbeam degradation below 3dB, from 
the figure, we should choose i?^ max /A < 1/2. This means that the maximum angle estimation error 
should satisfy 

*- £ 1r - 3P <79) 

and as R becomes large, the requirement of minimum angle ambiguity from (l79l becomes severe. 
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Fig. 14. Mainbeam degradation due to location estimation errors in the open-loop scenario. 



VII. Conclusion 

In this paper, we have analyzed the stochastic performance of random arrays for distributed collaborative 
beamforming, in the framework of wireless ad hoc sensor networks. It has been shown that under ideal 
channel and system assumptions, directivity of order N can be achieved asymptotically with N sensor 
nodes, as long as the sensor nodes are located sparsely enough. We have studied the average and the 
distribution of the beampattern as well as the distribution of the sidelobe peaks. Several forms of the 
CCDF of the beampattern have been derived and compared, with particular emphasis on the Gaussian 
approximation of the array factor. We have considered two scenarios of distributed beamforming and 
investigated the effects of phase ambiguity and location estimation error upon the resultant average 
beampatterns. 

Our main conclusion is that, given a number of nodes randomly distributed over a large disk, one 
may form a nice beampattern with narrow mainlobe and sidelobes as low as 1/N plus some margin 
for maximum sidelobe peaks. Also, the directivity approaches N if the nodes are located as sparsely as 
possible. However, our analysis is based on a number of ideal assumptions on the system and channel 
model. In practice, a number of open issues remain, such as applicability of beamforming when the 
destination or nodes in the cluster are in rapid motion or the channel suffers severe multipath fading. 
Also, specific algorithms should be developed for frequency offset correction of each node as well as 
methods for initial phase or location estimation. Finally, efficient protocols for sharing the transmit as 
well as calibration information among nodes are required. 



Appendix I 
Proof of Theorem[T] 



We first prove the following lemma: 
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Fig. 15. Function |2Ji(i)/i| 2 and its upper bound with x = 2.4445. 



Lemma 1: A generalized hypergeometric function 2F3 (I, § ; 1, 2, 3 ; — x 2 ) with i»l can be bounded 



as 



.an- 2 /,(i | ; 1,2,3; -^U I 



(80) 



where Co is a constant (co ~ 1.1727). 

Proof: We start with the integral form 



1 



Ji (x sin ■ 



x sin • 



de 



7T 







Vx 2 - t 2 



dt. 



(81) 



Since the asymptotic form of J\{x) given by dTol i is valid for x 3> 1, we have the following inequalities 

2 



Ji(t) 



t 



< cos (dot), for t < xq 



du < 



t+A 



— ^ du, for t > xq 



TTU° 



(82) 



(83) 



for some threshold value xq which should be determined numerically, and for some interval A > 0. The 
parameter a>o is chosen to be the smallest non-negative value such that 



cos(aoxo) 



irx 



3 ' 



(84) 



should hold, and this guarantees a continuity of the function at the threshold t = xq. Fig.1151 illustrates 
the relationship of d82l and d83l > with xo chosen as a cross point of t between the functions J\(t) and 
Y / 2/(7rt), yielding xq = 2.4445. The corresponding value of ao is 0.4664. Substituting (EU and d83l 
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into <f8TI) . we get for x > xq 



N 2 cos 2 (a t) 16 f x 1 
/ x < - / ~ T =^dt + - / ——=dt. (85) 

7T 7 VX 2 - t 2 7T^ A t 3 Vx 2 ~ t 2 



The first term on the RHS of dS3l is given by 

2 cos 2 (a t) 



^r{ i+ K^) 2+o(i/i4) } cos2( "° t>dt 

7tx \ 2ao / 



The second term on the RHS of ( I85t is given by 



i6 r i i6 



x 2xg 2x 3 } I V U/ I V^o 



4-2- f 1 " (-) 2 + (V* 4 )) + (M*)/* 3 ) 

7Hxq x \ V x / / v 

' * +0 (ln(x)/x 3 ) . (87) 



x vr^x, 



Consequently, we may write 



/(x) < I f xo + ^^o) + _LV +0 (ln(x)/x3) (88) 

and the second term on the RHS of (188b drops as x becomes large. With xo = 2.4445 and «o = 0.4664, 
the coefficient of 1/x can be calculated to be Co = 1.1727. ■ 
Proof: [Proof of Theorem Q From d3lT l. d32l . and Lemma [Q we have 

Ee. > £^L > \ = \ (Q.q\ 

For large N, the RHS of (g9j converges to O with /i = f£ » 0.09332. ■ 

Appendix II 

The Mean Number of Upward Level Crossings of a Gaussian Process 

In this appendix, we obtain the mean number of upward crossings of a given level of the zero mean 
Gaussian process based on the approach of Rice [10, 11]. Assume that X and Y are uncorrected zero- 
mean Gaussian processes with variance a 2 = a 2 = 1/2. Let u = sin (j^J and X' and Y' denote the 
corresponding processes differentiated by u. By assumption, X' and Y' become zero mean Gaussian 
processes. In order to calculate the variance, first consider the autocorrelation function of X at instants 
u = tii and v>2 given by 

px(ui,U2) = E z |cos (^zAnRu^j cos (zAitRu2^ I + other terms 

= -E t jcos (z4arR [u x + u 2 )) } + -E z jcos (zAitR ( Ul - u 2 )) } , (90) 
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where the other terms become zero by the zero mean assumption. Also for the same reason, the first 
term of the RHS of d90l) may be approximated by zero. Therefore, letting v = u\ — U2, we obtain 

p x {v) ±E Z {cos (WW) } . (91) 

Differentiating the above with respect to v twice, setting v = 0, and carrying out the statistical average 
with respect to z, the variance of X' is given by [12] 

a 2 x , = -p'x(0) = 2tt 2 R 2 . (92) 
Likewise, one may obtain a 2 , = a 2 ,, and the joint pdf of X, X' , Y, Y' is given by 

• /v v v v ' (x ' y ' x ' y } = <2kFo*o* exp r^i — 2^r J • (93) 

On changing the random variables in the polar coordinates via X = £1 cos 0, Y = O sin and integrating 
out © and 0', we obtain 

fapfau/) = we^-^e-^. (94) 



The number of positive (upward) crossings of the process uj at level a per interval du is given by [10, 
11] 

v(a)du = du / uj' fn,n>(a,u')ckj' = du2^y^rRae a . (95) 
Consequently, the mean number of upward crossings for the interval S3dB is given by 



which results in d54l) . 



E{u{a)} = du u(a), (96) 

Ju=sin(^),0ecS 3l iB 
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